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ABSTRACT 

We discuss cosmological hydrodynamic simulations of galaxy formation performed with the 
new moving-mesh code AREPO, which promises higher accuracy compared with the tradi- 
tional SPH technique that has been widely employed for this problem. In this exploratory 
study, we deliberately limit the complexity of the physical processes followed by the code for 
ease of comparison with previous calculations, and include only cooling of gas with a primor- 
dial composition, heating by a spatially uniform UV background, and a simple sub-resolution 
model for regulating star formation in the dense interstellar medium. We use an identical set 
of physics in corresponding simulations carried out with the well-tested SPH code GADGET, 
adopting also the same high-resolution gravity solver. We are thus able to compare both sim- 
ulation sets on an object-by-object basis, allowing us to cleanly isolate the impact of different 
hydrodynamical methods on galaxy and halo properties. In accompanying papers, we focus 
on an analysis of the global baryonic statistics predicted by the simulation codes ( Vo gels-| 
berger et aL]|2Ql and complementar y idealized simulati ons that highlight the differences 
between the hydrodynamical schemes ( [Sijacki et al.||2011 |). Here we investigate their influ- 
ence on the baryonic properties of simulated galaxies and their surrounding haloes. We find 
that AREPO leads to significantly higher star formation rates for galaxies in massive haloes 
and to more extended gaseous disks in galaxies, which also feature a thinner and smoother 
morphology than their GADGET counterparts. Consequently, galaxies formed in AREPO have 
larger sizes and higher specific angular momentum than their SPH correspondents. Interest- 
ingly, the more efficient cooling flows in AREPO yield higher densities and lower entropies 
in halo centers compared to GADGET, whereas the opposite trend is found in halo outskirts. 
The cooling differences leading to higher star formation rates of massive galaxies in AREPO 
also slightly increase the baryon content within the virial radius of massive haloes. We show 
that these differences persist as a function of numerical resolution. While both codes agree to 
acceptable accuracy on a number of baryonic properties of cosmic structures, our results thus 
clearly demonstrate that galaxy formation simulations greatly benefit from the use of more 
accurate hydrodynamical techniques such as AREPO and call into question the reliability of 
galaxy formation studies in a cosmological context using traditional standard formulations of 
SPH, such as the one implemented in GADGET. Our new moving-mesh simulations demon- 
strate that a population of extended gaseous disks of galaxies in large volume cosmological 
simulations can be formed even without energetic feedback in the form of galactic winds, 
although such outflows appear required to obtain realistic stellar masses. 

Key words: cosmology: dark matter - methods: numerical - galaxies: evolution - galaxies: 
formation - galaxies: haloes 



1 INTRODUCTION 

Cosmological simulations of galaxy formation provide a powerful 
technique to calculate the non-linear evolution of structure forma- 
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tion. In principle, such simulations can properly model the evolu- 
tion and formation of galaxies from first principles, but the enor- 
mous dynamic range and the many poorly understood aspects of 
the baryonic physics of star formation make this task extremely 
challenging. At present, the limitations imposed by the coarse nu- 
merical resolution achievable in practice and by the approximate 
treatment of the physics introduce significant uncertainties in sim- 
ulation predictions. Regardless, previous simulations have already 
proven instrumental for developing our current understanding of 



the formation and evolution of galaxies (e.g. Katz etal.|1992 jWein 



berg et al .'199T , Pearce et al.'1999','Kravtsov et al. 2002, S prmgel 
&Hernquist 2003a b, Borgani et al. 2004, Grain et al. 20091 



When only dark matter (DM) is considered, the newest gener- 
ation of cosmological simulation codes yields a consensus picture 
of important key results about the matter distribution in the Uni- 
verse, such as the large-scale distribution and the detailed internal 
properties of dar k matter haloes (Springel et al.^2005 . 2008^ Die- ^ 



mand et al.|2008||Boylan-Kolchin et al.|2009[ [Navarro et al.|2010| 



Klypin et al.| 2010V The computational methods employed in the 
codes vary between different implementations; however, the dis- 
cretization of the equations needed to follow the evolution of DM 
involves only gravity, and can only be done efficiently in terms of 
the N-body technique. Hence the requirements for obtaining con- 
verged and consistent results with different simulation codes are 
clear; all that is needed is an accurate gravity solver, a sufficiently 
large number of simulation particles and an accurate time integra- 
tion scheme. 

Over the last decades, a number of hydrodynamic codes have 
been developed and used in the field of computational cosmology 
(e.g.|Ryu et al.| 1993 [[Bryan et al.|1995l|Katz et al.|1996 Springel 
let al.[|2001| [Teyssier[[2002[ ']Wadsley et al. 2004, Springel 2005jT' 
They all account for the baryons in the Universe and try to follow 
their evolution self-consistently in a ACDM dominated universe by 
solving the equations of hydrodynamics coupled to gravity. How- 
ever, the discretization schemes and numerical algorithms used for 
this task differ widely from code to code, and unlike in pure DM 
simulations, large differences in the simulation outcomes are often 
found. Such simulations can be more directly related to observa- 
tions, unlike the pure dark matter simulations. 

In fact, hydrodynamic simulations have already explained sev- 
eral features of the observed large scale distribution of baryons, 
such as the existence and the properties of the Ly-a forest in the 
spectra of distant QSOs (e.g. Hernquist et al. " 1 996[ |M iralda-Escude [ 
|et al.[[7996| . In this, low density, low temperature regime differ- 
ent hydrodynamic techniques show good agreement between SPH 
and grid-based codes ( [Regan et al.[[2007| ). Numerical simulations 
played a crucial role in both our understanding of the properties of 
the absorbing gas and in the quantitative measurements of the mat- 
ter distribution in the Universe based on Lyman-a forest data (e.g. 
jCroft et al.|19"98l[Seljak et al.|2005[ ). 

Unfortunately, at higher densities characteristic of galaxies 
and their haloes, qualitative differences exist already at the simplest 
modeling level where radiative gas cooling and energetic feedback 
processes due to star formation are igno red (e.g.^ renk et al.ll999[ 
O'Shea et al.|2005[|Agertz et al.|2007||Tasker et al.|2008| fvazzal 



et al. 201 1| ). Such differences are sometimes apparent even be- 
tween codes using quite similar discretization techniques ( Springel[ 
|& Hernquist|2002 , Keres,2007j ). This emphasizes the need to find 
numerical methods that give the most accurate and reliable results. 
Understanding the origin and size of errors in a given computa- 
tional technique is clearly vital for fully exploiting the predictive 
power of cosmological simulations of galaxy formation. 



The most widely used methods to follow the evolution of 
baryons in cosmological simulations are Lagrangian smoothed par- 
ticle hydrodynamics (SPH) and Eulerian mesh-based hydrodynam- 
ics with adaptive mesh refinement (AMR). SPH has the great ad- 
vantage of being naturally adaptive, making it suitable for simu- 
lations with a large dynamic range in density, such as cosmolog- 
ical simulations where SPH automatically increases its resolution 
in collapsing regions. In addition, SPH can be easily combined 
with very accurate gravity solvers, which are often based on grav- 
itational tree algorithms (e. g. [Hernquist & Katz|1989| ) or particle- 
mesh-htree (e.g. Bagla'2002','Springel'2005') methods. SPH has ex- 
cellent conservation properties. In particular, it still manifestly con- 
serves total energy even when coupled to self-gravity, something 
not available in hydrodynamical mesh codes ( Miiller & Steinmetz| 
|1995| ). On the other hand, AMR codes offer higher accuracy in rep- 
resenting shocks, large density gradients and hydrodynamic fluid 
instabilities (e.g. [Agertz et al.|2007 ). 

It has become clear over recent years that both SPH and AMR 
codes suffer from weaknesses which are likely the primary cause 
for the disagreements found in some cosmological results obtained 
with these methods (e.g. Mitchell et al. 2009). For example, it has 
been shown that SPH in its standard implementation has limited 
ability to accurately follow fluid instabilities ( [Agertz et al.|2007| ). 
Furthermore, in SPH, shock capturing is done through the addi- 
tion of an artificial viscosity, which increases the dissipativeness 
of the scheme and can have unwanted consequences such as an- 
gular momentum transport in disks. Also, while the introduction 
of an artificial viscosity results in the correct properties of post- 
shock flow, the width of a shock in SPH can be significantly broad- 
ened (e.g. Hernquist & K atz 1989), potentially affecting the evo- 
lution of gas in complex cosmological environments (Hutchings 
& Thomas '2000^ for example causing in- shock cooling (^Creasey 



et al.|2011). On the other hand, Eulerian AMR codes are not free 



of unwanted numerical effects either. They may suffer from ad- 
vection errors, mesh-alignment effects, and potentially over-mixing 
( [Wadsley et al.|2008] ). In addition, AMR codes are intrinsically not 
Galilean-invariant, i.e. their truncation error depends on the bulk 
velocity relative to the grid. This can be a significant issue when 
the velocities of the flow are large and the resolution is limited (e.g. 
iSpringel 2010a , Robertson et al.|2010l l. 

Recently, [Springel[ ( |2010a| ) (SIO hereafter) proposed a new 
"moving-mesh" technique that promises accuracy advantages over 
both of these traditional approaches. This new method has been 
realized in the AREPO code and can be viewed as a hybrid of La- 
grangian SPH and Eulerian AMR, combining positive features of 
both methods while avoiding their most important weaknesses. In- 
stead of a stationary structured mesh as in AMR codes, AREPO 
uses an unstructured mesh given by the Voronoi tessellation of a 
set of mesh-generating points distributed within the computational 
domain. These points are allowed to move with the fluid, causing 
the mesh to be moved along as well, so that a natural adaptivity is 
achieved, just as in SPH. On the other hand, AREPO employs the 
same concepts of spatial reconstruction and Riemann- solver based 
flux computation as ordinary AMR codes, thereby retaining their 
high accuracy for shocks and fluid instabilities, and their low level 
of numerical dissipation. In addition, because the mesh can move 
with the fluid, the numerical solutions of AREPO are fully Galilean- 
invariant and have lower advection errors compared to codes with 
a stationary mesh. 

The AREPO code has already been successfully used in prob- 
lems of first star formation ( [Greif et "aLl|2011| ), and extensions 
for radiative transfer ( jPetkova & SpringeTl|2011| ) and magnetohy- 
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drodynamics ( [Pakmor et al.||2011| ) exist. Here we apply this new 
method for the fist time to cosmological galaxy formation sim- 
ulations in representative parts of the Universe. In what follows, 
we are primarily interested in clarifying to which extent the new 
hydrodynamic solver affects the properties of galaxies when com- 
pared to the widely employed SPH code GADGET (last described in 
Springel 2005 ). GADGET shares an identical high-resolution Tree- 
PM gravity solver with AREPO, allowing us to study differences 
caused by the hydrodynamics alone, unaffected by possible sys- 
tematics in the gravitational dynamics of the collisionless compo- 
nents, which were present in previous comparisons of mesh-based 
and SPH techniques (e.g. |0'Shea etal.|2 005 Hei tmann et al.|2 008|). 
Furthermore, we include an identical implementation of the basic 
physical processes that govern the evolution of the gaseous com- 
ponent in galaxy formation, such as radiative cooling and heating, 
star formation and feedback. Finally, we let AREPO carry out refine- 
ment and de-refinement operations on cells (if needed despite the 
moving mesh which normally keeps the mass per cell constant to 
good accuracy) in order to guarantee that the masses of cells never 
deviate significantly from a target given by the particle mass in the 
corresponding SPH comparison run. This ensures that we compare 
runs carried out with the two codes at about the same mass resolu- 
tion. 

In the last few years significant strides have been made in 
modeling disk galaxy formation with cosmological zoom-in sim- 
ulations based on the SPH technique (e.g. Gover nato et al.|20 "T0| 
[Guedes et al. 2011). These simulations often contain millions of 
resolution elements per halo and incorporate feedback processes 
that prevent an excessive collapse of baryons into galactic compo- 
nents. However, these successes in reproducing properties of ob- 
served galaxies have been reached for individual halos over a lim- 
ited mass range. It is not clear yet if the same model would be 
able to reproduce the observed population mix of galaxies over a 
wide range of properties. We stress, however, that our goal here 
is not to model galaxies with the highest possible resolution in 
a single system or to produce the most realistic looking galactic 
disks. Our comparison strategy is designed to detect systematic dif- 
ferences (that are projections of inaccuracies in used simulations 
techniques) for a large population of galaxies. While strong galac- 
tic feedback is an important component in galaxy formation and 
can induce large variations in the results (e.g. S cannapieco et al.| 
|2012| ), we here prefer an extremely simple subresolution treatment 
of feedback which can be implemented numerically in a well-posed 
and identical way, both in our moving-mesh code and in SPH. Our 
strategy attempts to ensure a clean comparison by avoiding a mod- 
eling of strong feedback processes, which would significantly com- 
plicate the interpretation of differences induced by the underlying 
hydrodynamic technique because parameterizing such feedback in 
a code agnostic fashion is quite difficult. 

In this Paper II of a series, we focus on an analysis of the prop- 
erties of galaxies and their haloes throughout cosmic time, and dis- 
cuss differences that the two simulation techniques imprint on the 
resulting population of galaxies. In |Vogelsberger et aTlpOl 1| here- 
after Paper I), we discuss global properties of baryons in different 
phases, the cosmic star formation history and various code charac- 
teristics, including technical details of our updated versions of the 
moving-mesh code AREPO and our SPH code GADGET. Finally, in 
a further companion paper ( [Sijacki et al.|201l| hereafter Paper III) 
we use idealized setups of situations relevant for galaxy formation, 
such as gas cooling in haloes and stripping of infalling galaxies to 
demonstrate where and how differences between the two codes oc- 



cur, which helps us to understand and interpret the differences we 
find in a fully cosmological environment. 

We note that the conclusions reached in this paper directly 
apply only to standard formulations of SPH (as described, for ex- 
ample in [Springel & Hernquist|[2002| . It is now well established 
that noisy and inaccurate pressure gradient estimates (for example 
see Read et al.|2010|), an absence of mixing at the particle level, and 
comparatively large artificial viscosity lead to a number of accuracy 
problems in SPH, in particular a limited ability to model fluid insta- 
bilities ( |Agertz et al.| 2007, Abel 20l T| and pr oblems in represent- 
ing subsonic turbulence (Bauer & Springel 2012 ), which are likely 
the most important factors contributing to the differences we iden- 
tify in this paper. We emphasise that this 'standard SPH' method 
or quite similar variants thereof have been used in the majority of 
the works on cosmological simulations of galaxy formation to date 
(both using large boxes and zoom-ins). A number of recent efforts 
proposed extensions or modifications of the traditional SPH formu- 
lation that aim to improve the accuracy of SPH (e.g. Wadsley et al. 
"Sons', "Price'2008';HeB & Springer2010;"Read et al. 2010 „Read& 
Hayfield 2012 , Abel 201 1 , Saitoh & Makino 2012 ). Our results do 
not necessarily apply to these new modified methods. They have 
yet to be tested for cosmological simulations of galaxy formation, 
and it remains to be seen whether they resolve the differences we 
find here in our comparison of standard SPH as implemented in 
GADGETS with a more accurate hydrodynamic scheme such as 
the one implemented in AREPO. 

The structure of this paper is as follows. In Section |2] we de- 
scribe our simulation setup, the physical processes included, the 
initial conditions, and the most important code modifications and 
methods used to identify haloes and galaxies. In Section [3] we 
present properties of galaxies and haloes in our moving-mesh and 
SPH simulations. In particular, we discuss the baryon content of 
haloes, galaxy mass functions, gas fraction and star formation rates 
of haloes, disk sizes and specific angular momenta of galactic 
baryons. We interpret and discuss our new results from the moving- 
mesh cosmological hydrodynamics in Section|4] 

2 METHODS 

In this paper, we describe only the essential setup information for 
our simulations, in the interest of brevity. We refer readers to the 
original code papers for a detailed discussion of GADGET ( |Springel| 
"2005 ) and AREPO (SIO), as well as to Paper I for technical infor- 
mation about specific parameters adopted for the codes and about 
their relative CPU- time performance. 

2.1 Simulation suite and initial conditions 

We construct the initial conditions at redshift z = 99 and evolve 
them to 2; = 0, adopting the following ACDM cosmology: Qrn = 
0.27, Qa = 0.73, Qb = 0.045, as = 0.8, = 0.95 and 
Ho = 70kms~^Mpc~^ (h = 0.7). These parameters are con- 
sistent with recent WMAP-7 measurements ( [Komatsu et al.|20TT| ) 
at the la level, except for Us which is consistent at 2 cr level. The 
initial power spectrum was approximated with the fit of |Eisenstein| 
|&Hul ( [T999l ). 

Our simulation suite follows a periodic box of size 
20 h'^Mpc at three different resolution levels: 2 x 128^, 2 x 256^ 
and 2 x 512^ particles/cells, respectively. Initially we start with an 
equal number of gaseous fluid elements and DM particles in both 
codes. While the number of DM particles stays fixed, the number 
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of gaseous elements changes over time. In GADGET this change 
owes to transformations of gas particles into star particles during 
star formation events, keeping the total number of baryonic resolu- 
tion elements constant. Likewise, in AREPO gaseous cells turn into 
star particles, but occasionally can also be refined and de-refined 
if needed to keep the mass resolution always close to the nomi- 
nal resolution present in the initial conditions. This leads to further 
changes in the number of cells, but overall, the total number of 
baryonic resolution elements stays approximately constant as well 
(see Paper I for more details). 

We note that we have used the same random phases and am- 
plitudes for overlapping modes on large scales for the different res- 
olution levels, so that our simulation set can be used for a resolution 
study of individual objects formed in our ACDM realization. In Ta- 
ble [T] we list the symbolic names of the simulation runs, together 
with their principal numerical parameters, such as the initial par- 
ticle/cell masses, the number of computational elements for all of 
our simulations, and the gravitational softening used for collision- 
less particles. 

2.2 Gravity, cooling and star formation 

The calculation of gravitational forces in all of our simulations is 
done with the algorithms described in |Springel| ( |2005| ). In short, 
we use a combination of a particle-mesh approach (e.g. |Hockney| 
|& Eastwood|1981i ) to calculate large-scale forces with the help of 
fast Fourier transforms, and a hierarchical multipole approximation 
(a tree algorithm, [Barnes &^ Hut 1986, I Hernquis t 1987 ) for short- 
range distances. This combination ensures a high computational 
speed nearly independent of the clustering state, and a uniformly 
high force resolution throughout the simulation volume. 

For the long-range PM computation, we employ a mesh 
size twice as large as the corresponding particle grid, meaning 
that for our L20nl28, L20n256 and L20n512 simulations, the 
mesh has 256^, 512^ and 1024^ elements, respectively. It is im- 
portant to emphasize that the gravitational force calculation is 
done in a consistent way in AREPO and GADGET, which enables 
us to cleanly identify differences caused by the hydrodynamics 
alone. The gravitational forces are softened at small distances as 
described in [Hernquist & Katz] ( |1989 ), with the corresponding 
Plummer-equivalent softening length fixed in comoving coordi- 
nates at 4/i~^kpc, 2/i~^kpc and l/i~^kpc, respectively. These 
values are applied to all collisionless particles in AREPO and all 
particle types in GADGET. Gaseous cells in AREPO use an adap- 
tive softening calculated as Cgas = 2.5 x (3yceii/47r)^/^ with an 
imposed minimum gravitational softening length equal to the fixed 
softening of our corresponding GADGET simulations. However, 
we have also run cosmological simulations where the gas cells used 
a fixed softening length as well, finding that this choice does not af- 
fect the properties of galaxies found in our runs in any significant 
way (see also Appendix of Paper I). 

We include all radiative cooling processes important for a hy- 
drogen and helium gas of primordial composition; i.e. line cooling, 
free-free emission, and inverse Compton cooling off the cosmic mi- 
crowave background. The rate equations are solved assuming col- 
lisional ionization equilibrium, as in [Katz et al.| { [1996] ). We also 
include a uniform ionizing background radiation field with the nor- 
malization and time dependence of jFaucher-Giguere et al.[ ( [2b09| ) 
that heats and ionizes the gas in an optically thin approximation 
( [Katz etal.|19 96^. 

Dense star-forming gas is treated via a simple two-phase sub- 
resolution model as in jSpringel & Hernquist[ ( [2003a| , which gives 



rise to an effective equation of state for the dense gas. We use 
the same parameters for the sub-resolution model as in| Springel &[ 
[Hernquist| p003aj ), which give a threshold density for star forma- 
tion of = 0.13 cm~^. The star formation timescale scales with 
density as oc and we set it to 2.1 Gyr at the threshold 

density in order to match the normalization of the local relation be- 
tween star formation rates and gas surface density of galaxies ( Ken-[ 
[nicutt| 1998 ). The conversion of gas mass to collisionless stellar par- 
ticles proceeds stochastically based on the estimated star formation 
rate. We simplify the original procedure by jSpringel & Hernquist[ 
(2003a) further so that each gas particle spawns only one genera- 
tion of stars. In addition, we modify the original model by letting 
gas above the star formation density threshold with a temperature 
higher than the effective two-phase temperature to radiatively cool 
in the ordinary way onto the equation of state, instead of applying 
the relaxation timescale of the multi-phase medium. 

We note that the pressurization of the star forming medium 
does not drive any outflow from the galaxy. For reasons of sim- 
plicity, we do not include additional violent feedback mechanisms 
(e.g. galactic winds) in the present study that could significantly 
lower the mass of a galaxy. Our model for supernova feedback can 
hence be viewed as a 'minimal feedback model' which cannot be 
expected to yield a particularly realistic galaxy luminosity func- 
tion, especially at the faint-end. In order to compare galaxy masses 
and luminosities of our predicted galaxy population with observa- 
tions we will consider simulations with galactic winds in AREPO 
in forthcoming work. We work here with this simplified model in 
order to allow us to identify differences between the hydro solvers 
as cleanly as possible. 

2.3 Cell regularization and mass homogenization in AREPO 

An important aspect for the accuracy and performance of AREPO is 
the use of a mesh-regularization scheme that ensures a reasonable 
quality of the Voronoi mesh, in the sense that unfavourable cells 
with very high aspect ratios are avoided. As described in detail in 
Paper I, we have adopted a new method where cells are only regu- 
larized if the maximum angle under which a face is seen from the 
mesh-generating point of a cell becomes large. Regularization here 
means that a small velocity component towards the geometric cen- 
ter of a cell is added to the motion of the mesh-generating point, 
which is otherwise just given by the gas velocity in the cell. 

In addition to the mesh regularization, we also employ on-the- 
fly mesh refinement and de-refinement operations, where appropri- 
ate, in order to ensure a constant mass resolution. In principle, this 
is not strictly required in AREPO, as the pseudo-Lagrangian nature 
of the scheme moves the mesh along with the mass, bringing more 
resolution to regions where this is needed most and maintaining 
roughly constant mass in cells. However, the resulting spread in 
mass tends to grow with time, and without explicit refinement and 
de-refinement measures, a relatively wide spectrum of cell masses 
may develop. An unwanted consequence would also be that star 
particles with considerable variation in mass would be created, 
which might lead to undesirable two-body heating effects. There- 
fore, we enforce cells to have a mass close to our desired target 
mass by splitting and merging individual cells if needed. In prac- 
tice, we allow cell masses to deviate by at most a factor of 2 
from the desired target mass. The latter is set equal to the SPH 
gas particle mass in the run with the corresponding initial resolu- 
tion, as listed in Table [T] so that we can ensure an equivalent mass 
resolution in the AREPO simulations. We note however that differ- 
ent refinement criteria other than mass could be chosen in AREPO, 
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Name 


Code 


Boxsize [{h~^Mpc)^] 


hydro elements 


DM particles 


^target/SPH[^ ^M©] 


mr>Mlh ^Mq] 


e [h ^kpc] 


A_L20n512 


AREPO 


20^ 


512^ 


512^ 


7.444 X 10^ 


3.722 X 10^ 


1 


G_L20n512 


GADGET 


20^ 


512^ 


512^ 


7.444 X 10^ 


3.722 X 10^ 


1 


A_L20n256 


AREPO 


20^ 


256^ 


256^ 


5.955 X 10^ 


2.977 X 10'^ 


2 


GX20n256 


GADGET 


20^ 


256^ 


256^ 


5.955 X 10^ 


2.977 X 10'^ 


2 


A_L20nl28 


AREPO 


20^ 


128^ 


128^ 


4.764 X 10'^ 


2.382 X 10^ 


4 


G_L20nl28 


GADGET 


20^ 


128^ 


128^ 


4.764 X 10'^ 


2.382 X 10^ 


4 



Table 1. Basic parameters of our simulation set. All calculations were performed in a periodic box with side length of 20 h~^Mpc. The number of hydro- 
dynamical resolution elements (SPH particles or Voronoi cells, respectively) becomes smaller with time due to star formation. In the AREPO simulations, 
additional re- and de-refinement operations are invoked to keep the cell masses close to a target gas mass mtarget, which we set equal to the SPH particle 
mass of the corresponding GADGET run. The comoving Plummer-equivalent gravitational softening length e is constant in GADGET, but adaptive for gas 
cells in AREPO. 



which has a similar flexibility as AMR codes in this respect. Fur- 
ther details on the technical implementation of refinement and de- 
refinement in AREPO are provided in SIO and Paper 1. 



2.4 Galaxy and halo identification 

We identify haloes in our simulations with the standard FOF ap- 
proach (Davis et al. 1985 ) with a linking length of 0.2 times the 
mean interparticle separation. The FOF algorithm is applied only 
to the DM particles, whereas baryonic particles/cells are later as- 
signed to their nearest DM particle and included in the halo to 
which the corresponding DM particle belongs. Self -bound con- 
centrations of mass, termed 'subhaloes', are then identified within 
each FOF halo using the SUBFIND algorithm ( |Springel et aL|2001| ), 
which has recently been modified to also account correctly for 
baryons in subhaloes ( |Dolag et al.|2009| ). 

The SUBFIND approach uses an adaptively smoothed density 
field (estimated from adding up the individual contributions of all 
particle/cell types) to find isolated overdense regions by means of 
an excursion set technique. Each substructure candidate identified 
in this way is then subjected to a gravitational unbinding procedure, 
leaving only the self-bound material in the subhalo. This effectively 
produces a list of subhaloes for each main FOF-halo. The bound 
part of the smooth background FOF-halo (after all embedded sub- 
structures have been removed) is defined as the 'main subhalo'. For 
the centre of each FOF halo (selected as the minimum of the grav- 
itational potential of the main subhalo) we also determine virial 
mass estimates with the spherical overdensity algorithm. Follow- 
ing a common convention in the literature, we adopt Mhaio=M2oo 
as estimate of the halo mass, where M200 is the mass contained 
within a spherical region of radius R200 around the centre with 
overdensity of 200 with respect to the critical density. We use the 
term 'central galaxies' for the galaxies in the central parts of the 
main subhalo of a given FOF halo. 

Galaxy identification can be a complex task, especially in 
GADGET, where a large number of cold gaseous clouds form in 
massive haloes at large radii ( [Kaufmann et al.||2009[ |Keres et aT] 
|2009a[ IPuchwein et al.pOlO i SUBFIND tends to associate those 
with the central galaxies, which can skew some of their proper- 
ties, like their angular momentum, half-mass radius and, in extreme 
cases, even their total mass. In order to protect against this effect 
and to make the galaxy identification closer to what would be in- 
ferred observationally, we adopt a few further criteria after sub- 
halo identification in order to select the galactic baryons and con- 
struct the galaxy catalogue. Specifically, we assume that galaxies 



include only reasonably cold gas with temperature T < 30, 000 K, 
apart from the star-forming gas at high-density, whose effective 
temperature can be somewhat higher (due to the pressurization 
of the two-phase medium). We also require gas to have bary- 
onic overdensities higher than 2000 which allows us to include 
all of the star forming gas but also extended gaseous disk out- 
side of the star-forming region. Furthermore, we apply a radial cut 
of r < 30kpc/h comoving to both gas and stars to avoid that 
gaseous clouds and filaments in the outskirts of haloes, or sub- 
stantial numbers of intra-group/intra-cluster stars ( [Puchwein et aL| 
|2010| ) are included in our definition of a galaxy. While these cri- 
teria are clearly crude, they are motivated by observed properties 
of galaxies that contain mostly cold dense gas and stars. To ex- 
clude poorly resolved galaxies we concentrate on galaxies more 
massive than 256 x mtarget, or when selecting by halo mass, we 
use M200 = 256 X mtarget X Qrn/^b as our adopted resolution 
limit. 



3 RESULTS 

3.1 Baryon content of haloes 

In Paper I, we have shown that the dark matter halo mass functions 
agree very accurately between our different simulation techniques, 
and also converge well as a function of numerical resolution. 
Here we concentrate on the baryonic properties of these haloes. In 
Figure[T] we show the mass fraction of baryons and of gas in haloes 
as a function of parent halo virial mass, comparing AREPO with 
GADGET at different resolutions. Both fractions are expressed in 
units of the value /Qm for the universe as a whole, which would 
be the expected mass ratio if the baryons would perfectly trace 
dark matter. Figure [T] illustrates that this is only approximately 
true, because the baryonic fraction is a function of halo mass and 
often deviates significantly from unity. 

In fact, in low mass haloes below ^ 10^° /i^^M© at z = 0, 
or 5 X 10^ /i~^M0 at 2; = 2, the fraction of collapsed baryons 
is significantly lower than the universal value. This is a well known 
consequence of the IGM photo-heating by the UV background, ef- 
fectively weakly "pre-heating" the gas, modifying the peaks in the 
cooling curve ^Katz et al.||1996j ) and preventing collapse and ac- 
cretion of gas into low mass haloes (e.g. |Efstathiou|199^|Thoul &| 
|Weinberg|[T996l |Gnedin|[2000l ). For the UV background intensity 
used in this paper a rapid decline in baryonic fractions is in fact 
expected below - 5 x 10^ H'^Mq at z = 2 and - 10^° H'^Mq 
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Mhco [M© /h] M,„„ [Me /h] 

Figure 1. Baryonic content of haloes for our simulation set at z = 2 (left panels) and z = Q (right panels). The lines show the total amount of baryons in 
haloes normalized by the baryon to total matter ratio in our adopted cosmology. The panels on top compare the total baryonic content in haloes of AREPO 
simulations (shown in red) with those in corresponding GADGET runs (blue lines). Different resolutions are displayed with different line styles, as labeled. In 
the bottom panels, we show the same comparison for the baryons in gaseous form (i.e. excluding stars), normalized the same way. The different vertical lines 
indicate a mass of 256 x mtarget x Q^rn/^h for ^^ch resolution. 



at z = (see also |Hoeft et al.|2006[[O'kamoto et al.|2008[|Faucher-| 
IGiguere et al.|201l 1). The mass scale at which haloes transition into 
the baryon-poor regime is consistent between SPH and AREPO and 
slightly increases when the resolution is poor (see the 2 x 128"^ 
runs). 

More massive haloes have a baryonic content relatively close 
to the universal value, but with interesting trends as a function of 
halo mass, redshift and simulation technique. At z = 2, the bary- 
onic fraction peaks 15% above the universal value at halo masses 
around 2 x 10^° /i^^M© in both codes. However, while in GAD- 
GET this drops again to about unity for Mhaio > 10^^ H'^Mq, in 
AREPO there is still a baryon excess of around 10% in such mas- 
sive haloes. Similarly, at z = the peak in the baryonic fraction 



has moved up a bit in halo mass scale to ^ 10^^ /i^^M© and is 
even a bit higher. At this epoch, well-resolved AREPO haloes can 
have a 25% excess of baryons while in GADGET the peak is slightly 
lower. In more massive haloes, Mhaio > 10^^ /i^^M©, the total 
baryon content drops again to around unity in SPH but remains 
around 10% higher than the universal value in AREPO. This sys- 
tematic difference in massive haloes is very robust to changes in 
resolution. While this result may at first appear a bit puzzling, we 
note that such a systematic offset between the baryon fractions in 
SPH and mesh-based hydro codes has in fact been noticed earlier. 
Already the non-radiative simulations of the Santa Barbara cluster 
comparison project (Frenk et al. 1999 ) found that SPH runs had 
only ^ 90% of the expected baryons, whereas grid codes found a 
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ratio close to unity (see also |0'Shea et al.|2005||Crain et al.|2007| ). 
If radiative cooling is included, simulations have been found to pro- 
duce a higher baryon fraction at any given radius, including out to 
the virial radius (e.gJKravt sov et al.|20()5l|Ettori et al.|2006j ). Our 
results are thus consistent with these earlier studies, even though 
we find a somewhat larger difference between SPH and mesh-codes 
compared to Kravtsov et al. ( 2005 ) when cooling is included. Note 
however that our small simulation box does not probe the mass 
scale of galaxy clusters examined in |Kravtsov et al.|f2005| ), where 
the effect is likely smaller. 

We speculate that the origin of the difference between SPH 
and the mesh codes ultimately lies in the higher dissipative heat- 
ing rate we found for SPH haloes in Paper I, which effectively 
leads to a higher pressurization of the baryons in the gaseous at- 
mospheres, thereby pushing them out slightly in SPH relative to 
the mesh-codes. When cooling is included, the baryons removed in 
the halo centres allow the outer gas to slide in, thereby increasing 
the baryon fraction inside the virial radius. Since the large AREPO 
haloes cool more gas, this effect is stronger in the mesh code. Also, 
the reduction in strength of this effect towards higher masses can be 
understood as a consequence of a declining cooling efficiency. In 
fact, when our results are extrapolated to cluster scales, they appear 
entirely consistent with findings in galaxy cluster simulations that 
include radiative cooling, and in particular, with the weak trend of 
a falling baryon fraction with increasing mass found there ( Borgani 
|& Kravtsov||2009l ). Before moving on we remark that the differ- 
ence in baryon fraction also means that there is a slight but sys- 
tematic difference in the virial masses of matching large haloes in 
AREPO simulations compared to GADGET. This difference is how- 
ever < 10% in most cases which means that the differences between 
the codes presented as bin averaged results at a given halo mass are 
essentially the same as they would be for binned matched samples 
of halos. 

In the bottom panels of Figure [T] we consider the gas frac- 
tions of haloes, at redshifts z = 2 and z = 0. The gas fraction of 
haloes is quite low, amounting to ^ 40% of the universal baryon 
fraction for large haloes, implying that most of the baryons in these 
haloes are actually locked up in stars. In fact, in our high resolu- 
tion simulations at 2; = this is true for all haloes above the UV 
background heating limit. This extremely efficient conversion of 
gas into stars is a consequence of the absence of very strong feed- 
back processes in our simulations that could eject material from 
the bottom of the potential well of haloes in the form of winds or 
superwinds. Instead, the amount of stars formed is primarily reg- 
ulated by the cooling efficiency of haloes, and the differences we 
find between the simulation techniques reflect differences in this 
cooling efficiency. In particular, there is less gas in large haloes at 
late times in AREPO, despite the fact that these haloes contain more 
baryons over-all compared with GADGET. This is consistent with 
the stronger cooling flows and higher star formation rates we found 
in these objects (see also Paper I). On the other hand, in small mass 
systems just above the cooling threshold, AREPO retains a bit more 
gas, both at 2; = and 2; = 2, suggesting that here the mesh code is 
forming stars somewhat less efficiently. In this regime, the moving- 
mesh code is also more sensitive to numerical resolution than SPH, 
showing a stronger increase in the remaining gas fraction if the res- 
olution is degraded. As we will discuss later on in the paper, this 
is in part a consequence of differences in disk sizes of galaxies and 
the inability of low resolution simulations to resolve high density 
gas, both affecting the star formation efficiencies. 



3.2 Galaxy mass functions 

Figures [2] and |3] show galactic baryonic mass functions (stars plus 
gas associated with galaxies) and galactic stellar mass functions 
at 2; = 2 and 2; = for all of our simulations. For our highest 
resolution simulations, the galaxy mass functions are in relatively 
good agreement at 2; = 2 for both simulation techniques, up to 
masses of several 10^° /i^^M©. At the high-mass end, the masses 
of AREPO galaxies are however systematically higher than those of 
GADGET galaxies, by a factor of ^ 1.5 - 2. 

The agreement at the low mass end, is arguably to be expected. 
Here the galaxies accrete most of their gas rapidly with rates com- 
parable to the gas accretion rates onto haloes ( jKeres et al.||2005[ 
|Ocvirk et al.||2008] ). Furthermore, in this low halo-mass regime, 
smooth accretion dominates over growth by mergers with other 
galaxies ( |MuraH et al.p002l|Keres et al.|2005|[Guo & White|20Q8| l, 
owing to rapid gas infall and a limited mass range available for pro- 
genitors close to the UV-background suppression mass. Therefore, 
there is comparatively little opportunity in these low mass haloes 
that differences in the treatment of hydrodynamic interactions sig- 
nificantly affect the growth of galaxies. However, we note that this 
is only approximately true as there are clearly some differences 
when then resolution is low, as we discuss below. 

At 2; = 0, some more prominent differences become appar- 
ent. While the highest resolution simulations again show reason- 
ably good agreement in the galaxy mass function below 2 — 3 x 
10^° /i-^M0,at higher masses, the differences in the galaxy mass 
reach a factor of ^ 2. While both numerical methods converge to 
the same galaxy mass function at the low mass end as the resolu- 
tion is improved, at the high mass end, AREPO and GADGET still 
converge individually, but they do so to a different converged result. 
This systematic difference manifests itself in higher galaxy masses 
in AREPO, which we argue is linked to a stronger cooling efficiency 
in massive haloes, an issue we will return to in Section |T6l 

We find that GADGET simulations show convergence in the 
galaxy mass function for galaxies more massive than ^ 250 bary- 
onic particles (about 8 x A^ngb for our runs). To compare this res- 
olution limit to the convergence achieved with AREPO, we use ver- 
tical lines in Figures [2] and |3] to mark the mass scale that corre- 
sponds to 256 times the target cell mass. AREPO appears to con- 
verge more slowly, requiring up to ^ 5 times higher mass at 2; = 2. 
At 2; = 0, the resolution dependence of both codes is even stronger 
and the convergence of the high-resolution simulations is somewhat 
poorer than at 2; = 2, likely caused by differences in resolving the 
halo scale affected by the UV background. We caution that in the 
cosmological simulations without galactic outflows improved res- 
olution always leads to a higher fraction of baryons locked up in 
stars, especially before reionization, as there is no efficient process 
to regulate formation of lower mass objects. Furthermore details 
of galactic structure will depend, in part, on resolution, especially 
since the hierarchical buildup proceeds from smaller objects that 
are always less well-resolved. Therefore, precise criteria for con- 
vergence depend on the galactic property one is interested in. Note 
also that faster convergence of SPH for some of the properties does 
not guarantee that it is a more accurate numerical scheme. In fact, 
in principle it is possible that one may even converge to a wrong 
result (e.g. Springel 2010b ). This is because in practice once typi- 
cally only increases the total particle number and spatial resolution 
in an SPH convergence study. However, an increase in the number 
of neighbours would also be needed for convergence in order to 
progressively beat down the SPH noise (Rasio 2000^ ), as we discuss 
in detail in Paper I. 
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Figure 2. Baryonic mass function (gas+stars) of galaxies identified in our simulation set di z = 2 and z = 0. Runs with AREPO are shown as red lines, those 
with GADGET as blue lines. The different resolution levels are indicated with different line styles, as labeled. The vertical lines indicate the mass of 256 x 
target cell/particle mass. At high redshift and at high resolution, the mass functions are similar between the codes. At z = 0, the differences are larger, as we 
discuss in the text. 
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Figure 3. Stellar mass functions of galaxies for our whole set of simulations at 2; = 2 and 2; = 0. As in Figure [2] we use red lines for AREPO and blue 
lines for GADGET, with different line styles corresponding to different resolutions, as labeled. The vertical lines indicate a mass of 256 x mtarget for the 
corresponding resolution level. We note that the differences in the stellar mass function are slightly larger than the differences in the galaxy mass function, 
caused by systematic offsets in the gas fractions of galaxies in AREPO and GADGET simulations. 



The galaxy stellar mass functions shown in Figure [3] give a 
qualitatively similar picture, if anything, the differences between 
the codes are even more evident. Here the larger offsets between 
the codes at low masses are caused by differences in galactic gas 
fractions at a given galaxy mass, as we discuss more in Section [34| 
At the high mass end, where stars form the bulk of galactic mass, a 
strong systematic offset in stellar masses is again clearly evident. 

To summarize, while some moderate differences exist in low 
mass galaxies, these cannot be unambiguously attributed to sys- 
tematic effects of the hydrodynamical solvers, owing to numerical 
resolution and galaxy finder uncertainties in this regime. On the 
other hand, the differences we find between the moving-mesh and 



SPH codes at the massive end are very significant and robust. Here 
AREPO produces more massive galaxies, which is one of the key 
differences we will return to several times in the rest of the paper. 



3.3 Relation between galaxy and halo masses 

Some of the trends we inferred from the galaxy and stellar mass 
functions are more explicitly demonstrated in Figure |4] where we 
show the stellar mass and baryonic mass for central galaxies as a 
function of their parent halo mass M200. The measurements make 
it clear that the agreement of the baryonic masses of galaxies in our 
two simulation techniques is indeed excellent in haloes of mass 
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Figure 4. Dependence of galaxy and stellar mass of central galaxies on halo mass at 2; = 2 and z = 0. Upper panels show the galaxy mass while the lower 
panels show the stellar mass. We compare results for AREPO (red lines) and GADGET (blue Hnes) at different resolutions, as labeled. It is clear that at high 
resolution the galaxy masses of AREPO and GADGET galaxies are indeed very similar in haloes with M^aio ^ 10^^ h~^M(7), but in more massive haloes 
the galaxies in the moving-mesh simulations are systematically more massive. This relative mass difference increases slowly with halo mass. Vertical lines 
show 256mtarget ^m/^b at different resolutions. 



M200 < lO^^/i ^M0, whereas galaxies in more massive host 
haloes are systematically more massive in AREPO runs. 

Because our galaxies have a relatively low gas fraction at 
z = 0, the stellar and baryonic masses are comparable at all halo 
masses at late times. In contrast, at 2; = 2, a lower amplitude of the 
stellar mass function when compared with the galaxy mass func- 
tion is apparent, caused by a much higher gas content of galaxies 
at a given halo mass at this epoch. The differences in the galaxy 
masses at low resolution are more dramatic which is caused by the 
sharply declining masses of poorly resolved galaxies in low mass 
host haloes. In this case, simulations are not able to properly resolve 
galactic densities any more; in fact, here the gravitational softening 
can become comparable to the galaxy size, so that the collapse of 



baryons in low mass haloes is seriously affected. We note that our 
findings also confirm previous results that the amount of baryons 
accumulated in central parts of haloes substantially exceeds the ob- 
servationally inferred galaxy masses (e.g. |Guo et al.|2010}|Behroozi| 
|etalj2010) if no efficient feedback mechanism is invoked that ei- 
ther ejects gas from these regions or somehow prevents the cooling 
in the first place ( |Keres et al.|2009b| ). At the massive halo end, the 
AREPO simulations make this problem even more acute, strength- 
ening the arguments for the need of very strong feedback mecha- 
nisms. 



© 0000 RAS, MNRAS 000, 000-000 



10 D. Keres et al. 





Figure 5. Gas fractions of central galaxies, Mgas/Mgas+stars, as a function of halo mass for our simulation set (left panels). We use red lines for AREPO and 
blue lines for GADGET, with different line style corresponding to different resolutions, as labeled. The gas fractions are systematically higher in AREPO than 
in GADGET. We also show the fractions of galactic gas that reach star-forming densities MsFgas/^gas (right panels) to illustrate how this depends on halo 
mass, numerical technique and resolution. Even though both quantities are resolution dependent, the overall systematic trends clearly show higher galactic gas 
fractions and a lower relative fraction of gas in star-forming regions in AREPO. 



3.4 Galactic gas fractions 

In the left-hand panels of Figure [5] we show the gas fractions 
Mgas/ (Mgas + Mstars) of Central galaxies in our simulation runs 
at z = 2 and z = 0. Lines give median gas fractions at a given 
parent halo mass at the different resolutions. At any given resolu- 
tion, the gas fractions are higher in AREPO than in GADGET runs, 
but the magnitude of the difference depends on the resolution of 
the simulation. It is interesting that even in our highest resolution 



simulations, the gas fraction of AREPO galaxies are systematically 
higher than in GADGET indicating robust differences in galaxy evo- 
lution between the codes. 

The fact that we find some resolution sensitivity is not too 
surprising as in our lowest resolution runs the gravitational soften- 
ing can be a sizable fraction of the virial radius of low mass halo 
haloes and is comparable to the sizes of the corresponding galaxies. 
The higher resolution simulations are expected to better resolve the 
dense gas in galaxies, leading to a higher star formation rate, and 
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Figure 6. Projected gas distributions of two interacting galaxy systems in GADGET (left panels) and AREPO (right panels) at redshifts z = 2. The panels 
show regions of lOOkpc/h comoving on a side. Upper panels show projections of gas in a major merger of two low mass galaxies while lower panels show a 
more complex merger involving three galaxies. The tidal features, tails and bridges, seen in the gas are much more prominent in AREPO, demonstrating that 
the galaxies simulated with the mesh code contain a larger amount of cold, rotationally supported gas. 



consequently to a lower residual gas fraction at a given halo mass. 
Once the relevant densities are properly resolved, converged values 
for the gas fractions should be reached, reflecting a balance be- 
tween gas supply and the depletion of gas by star formation. How- 
ever, the properties of the galaxies and of their gas content will 
still be affected by their hierarchical formation path. The buildup 
of structure from quite poorly resolved small sub-components can 
therefore impact the gas infall rates, and also the density structure 
of much larger galaxies. It is hence not sufficient for convergence 
to achieve a sufficiently high resolution to adequately represent a 
given galaxy, one also needs to represent well enough a good part 
of the merger history leading to this object. This makes hierarchical 
galaxy formation a particularly difficult numerical problem. 

It is clear that finite spatial resolution and the limited ability 
to follow the full hierarchical buildup will affect both of the codes 
in some ways. One might however expect numerical problems rele- 
vant to the density structure of dense galactic gas to be more severe 
in SPH. This is because SPH relies on an artificial viscosity in order 
to capture shocks and to produce the associated entropy increase. 
It has been shown that the artificial viscosity can lead to unwanted 
angular momentum transfer in the gas of a rotating galactic disk, es- 
pecially when the disk is sampled with a small number of particles 
( [Navarro & Steinmetz|1997||Okamoto et al.|2003l ). The net effect 
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is that the disk mass looses angular momentum to material further 
out in the halo, making these systems more compact in GADGET, 
resulting in turn in a faster conversion of baryons into stars and in 
lower residual gas fractions. 

In Figure [6] projected gas densities for two different interact- 
ing galaxy systems in GADGET and AREPO are shown, supporting 
this picture. The tidal features, such as bridges and tails seen in 
the gas are prominent and well defined in AREPO, demonstrating 
that the galaxies simulated with the moving-mesh code contain a 
large amount of cold, rotationally supported gas ( [Toomre & Toomre| 
1972, D'Onghia et al. 2010). These features are much less promi- 
nent for the same systems in GADGET. This is also consistent with 
our findings in Section [3J] where we show that even in low mass 
haloes GADGET galaxies are more compact than AREPO galaxies. 
We note that since poorly resolved galaxies are also the progenitors 
of better resolved, more massive galaxies, this effect can propagate 
to larger systems and induce systematic differences in them. 

The gas fractions in galaxies decrease with time in all simu- 
lations, irrespective of the resolution and the code employed. The 
bulk of this trend is caused by a decreasing rate of gas supply into 
galaxies at a fixed halo mass, as commonly predicted by semi- 
analytic models and cos mological simulations (White & Frenk 
|T99T|purali et aL|2002||Hernquist & Springel|2003||Keres et" 
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|2005[|Dekel et al.|2009 l). Owing to the short (observationally mo- 
tivated) star formation timescale in our ISM sub-resolution model, 
the rate of conversion of dense gas into stars is practically con- 
trolled by the rate of gas infall. However, this is correct only to the 
extent that most of the galactic gas is at high enough densities to 
form stars, which depends on the disk sizes. Differences in disk 
sizes will therefore influence the effective star formation timescale 
and the gas fractions. The characteristic sizes of galactic disks at a 
fixed halo mass increase with time (see below), enabling a larger 
amount of gas to reside in a disk for a given star formation rate. 
This slows down the evolution of the gas mass in the disk, result- 
ing in a relatively modest drop in gas fractions with time (see also 
IDutton et al.TM^ . 

In the right-hand panels of Figure[5] we also show the fraction 
of galactic gas that resides in the star forming two-phase medium, 
i.e. Mgas,SF/Mgas, as a function of halo mass. At z = 2, in cen- 
tral galaxies of low mass haloes, there is more gas outside of the 
star forming region in AREPO than in GADGET. At z = 0, we find 
a similar difference, except for the most massive galaxies where 
the fraction of star-forming gas is comparable between the codes. 
Given that AREPO galaxies have also higher gas fractions, this 
implies that their galactic gas is more extended, allowing a large 
quantity of gas to stay below the star formation threshold. We con- 
firm this in Section [3J] where we demonstrate that galaxies in the 
moving-mesh code are more extended. It is interesting to note that 
in ~ 10^^ /i~^M0 haloes at z = the gas fractions in AREPO are 
a factor of ^ 2 higher than in GADGET, even though differences 
in the fraction of star-forming gas are relatively minor. As it turns 
out these galaxies are also accreting much more gas from their sur- 
rounding hot haloes (see Section |3.6| ), which contributes to their 
higher gas fractions. 



3.5 Star formation rates of simulated galaxies 

In Figure |7] we show the dependence of the star formation rate of 
central galaxies on their parent halo mass at redshifts z = 2 and 

z = for both codes and at three different resolutionO 

The most obvious trend in Figure |7]is that the SFR is a strong 
function of the parent halo mass. This dependence is stronger in the 
case of AREPO with almost linear dependence on halo mass for the 
10^° h~^MQ < Mh < 3 X 10^^ H'^Mq range. It is clear that 
for GADGET this relation is shallower, especially in haloes above 
3 X 10^^ h~^MQ where at z = the SFR is only weakly de- 
pended on halo mass. Given the close to linear relation of galaxy 
mass and halo mass in Figure |4] this also implies a steep relation 
between SFR and galaxy mass over an extended mass range, sim- 
ilar to observed trends (e.g. [NoeskeeT al. 2007, Dadd i et al.|2007[ 
[Salim et al.|2007| ). We find that in AREPO this relation has a higher 
normalization and stronger dependence on halo mass at the mas- 
sive end than in SPH simulations. This might help to explain the 
high observed star formation rates of galaxies at 2; = 2 — 3 ( |Forster| 
[Schreiber et al. |2009 ), which appear to form stars at rates higher 
than in previous models (Dave 2008 ). SFRs of simulated galaxies 
at a given mass show a rapid decrease with time in both codes, 
qualitatively consistent with observed trends and previous findings 



Note that the information displayed in this figure differs from Figure 10 
in Paper I by restricting the star formation measurement to the central 
galaxy in a halo instead of showing it summed for all galaxies within the 
virial radius of a halo. 



from cosmological simulations (e.g. |Keres et al.|2009a||Wuyts et aL| 
I20TT] ). 

At high redshift, GADGET galaxies show good convergence in 
the SFRs of large haloes, while AREPO exhibits a slightly stronger 
residual drift with resolution. However, for both methods, one typi- 
cally finds that low resolution causes a reduction of the SFR in low- 
mass objects, but leads to an overestimate of the SFRs in high-mass 
objects. This can be understood as a generic consequence of hier- 
archical galaxy formation. If limited resolution suppresses some 
of the star formation in low-mass progenitor objects, more gas is 
left over to support a higher SFR in large objects once they have 
formed. 

The SFRs in GADGET and AREPO reach similar values for 
low-mass parent haloes in our highest resolution simulations, 
which is a consequence of the relatively simple gas infall physics 
in this regime. However, it is clear that there is a substantial sys- 
tematic difference in massive haloes where AREPO shows a much 
higher star formation rate in the central galaxies. The difference 
is especially acute for the mass range Mhaio ~ 5 x 10^^ — 3 x 
lO^^/i~^M0 at z = 0, where the SFRs of central galaxies in 
AREPO are a factor of 3-5 times larger than for the same galax- 
ies in GADGET. This is also the reason why the global star forma- 
tion rate in AREPO is significantly higher at low redshift than in 
GADGET, as discussed in Paper I. Since the SFR in these massive 
haloes is largely regulated by the gas supply, and late-time mergers 
are relatively gas-poor, the main cause of these high SFRs at z = 
in AREPO must lie in a more efficient cooling in massive haloes. 
Similar differences are apparent for masses Mhaio ^ 10^^ Mq at 
z = 2, where AREPO 's SFRs are a factor of 2-3 higher than in 
GADGET simulations. This is a clear example of dramatic changes 
in the properties of galaxies when a more accurate hydrodynamical 
method is used to model galaxy formation. 

Interestingly, at high resolution, the SFRs of central galaxies 
in our most massive haloes, Mhaio ^ 10^^ /i^^M© are quite com- 
parable between the codes and are relatively low relative to their 
halo mass. In the next subsection we explore the hot gas profiles in 
massive haloes to examine these differences more explicitly. 

3.6 Gas profiles of massive haloes 

In Figures [8] and [9] we show average density, temperature, entropy 
and radial velocity profiles of-- lO^^/i~^M0 and-- 10^^ H'^Mq 
haloes at 2; = 0, respectively. To measure these profiles, we first se- 
lect all haloes in a narrow mass range as labeled in the temperature 
panels. For each halo, we find its spherically averaged profiles for 
non star-forming gas with T > 10^ K. We then stack the haloes by 
finding the median value over all systems in the given mass range. 
This is to avoid local distortions owing to infalling substructures 
that often have slightly different locations in AREPO and GADGET, 
and hence helps to focus attention on the smooth hot halo gas only. 
We note that all the haloes in this mass range are dominated by the 
hot halo gas. To allow a more faithful stacking, we express the ra- 
dius in units of R200, density in terms of enclosed baryonic density 
within R200, velocity in units of V200 = y^GM2oo/R2oo and tem- 
perature as a virial temperature, Tvir = O.S/imH'^^ioo/^B, where 
TTiH is the proton mass, /i is mean molecular weight of the gas, 
which we here take to be fully ionized, and /cb is Boltzmann con- 
stant. 

There are interesting systematic differences in the hot halo 
structure between the different simulation codes, that are especially 
prominent for our lower mass sample. The density profiles are 
steeper in AREPO compared with GADGET, showing at large radii 
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Figure 7. Star formation rates of central galaxies as a function of their parent halo mass for all of our simulations. Blue lines show GADGET galaxies while 
red lines show AREPO galaxies a.t z = 2 and z = O.ln low mass haloes in our highest resolution simulations galaxies form stars with similar rates. However 
in haloes more massive than 2 x lO-*^-*^ h~^MQ AREPO's star formation rates are higher at a given mass. These differences reach a factor of ~ 5 in Milky 
Way-mass haloes at 2; = 0. Lines show medium star formation rates at a given mass with different line types corresponding to different resolutions as labeled. 



lower gas densities, and in the central parts higher density. The tem- 
perature profiles are quite similar between the two simulation tech- 
niques for R > 0.2i?2oo. However, the inner region clearly shows 
a decrease in gas temperature in AREPO simulations indicative of a 
cooling flow, while the temperature continues to rise at small radii 
in GADGET. At lower resolution, there is an even stronger decrease 
in temperature and an increase in density in the innermost halo re- 
gion in AREPO. Even GADGET shows some mild decline of the cen- 
tral temperature and a density increase, albeit with a much lower 
magnitude and confined to smaller radii. These density and tem- 
perature differences are reflected in corresponding trends for the 
entropy of the halo gas. The entropy in moving-mesh simulations 
is significantly higher in the halo outskirts but is lower in the central 
~ 0.1 R200 of the halo. Finally, the mean radial streaming velocity 
of the gas in Figure|8]shows a larger negative amplitude for AREPO 
than for GADGET. This can be viewed as a tell-tale sign of a larger 
cooling rate out of the AREPO haloes. 

These trends provide a clear explanation for the large differ- 
ences in SFRs between the two simulation techniques visible in 
Figure [7] which are induced by a much more efficient cooling in 
haloes with mass ^ 10^^ /i^^M© in AREPO. In the more mas- 
sive, group- sized haloes shown in Figure[9] the differences between 
the gas profiles are however much smaller. There are still similar 
trends of a lower gas density in the outskirts and higher densities 
and lower entropies in the central parts of AREPO haloes that in- 
dicate some residual cooling difference. However, the temperature 
profiles are very similar between the codes, except for a slightly 
steeper rise of the central temperatures in the moving-mesh ap- 
proach, which is likely caused by a deeper potential well owing to 
a much higher accumulation of baryons in the centre. Also, in both 
codes the central entropy remains relatively high and the central 
temperature does not show signs of strong cooling. This explains 
the relative low SFRs over this halo mass range and the relative 
consistency between AREPO and GADGET. The long cooling times 
are also reflected in the radial infall profile, which shows that the 



hot halo gas is essentially in hydrostatic support and has very lit- 
tle average motion in the radial direction. We note that the high 
entropies in massive haloes in GADGET are consistent with pre- 
vious results from cosmological simulations with cooling and star 
formation ( [Keres et al.|2009a| ), but interestingly, the trend we find 
here between AREPO and GADGET is opposite to what is typically 
found in non-radiative simulations of cluster formation where the 
central entropies of massive haloes are found to be higher in mesh 
codes ( Frenketal..l9 99). 

Interestingly, we also find that the acceleration profiles due 
to hydrodynamical pressure forces (— VP/p) are very similar be- 
tween the two codes, as expected for the same haloes in quasi hy- 
drostatic equilibrium, because the pressure forces need to balance 
gravity, which is dominated by the same dark matter distribution in 
both codes. Is is then at first a bit puzzling to understand why the 
AREPO haloes are cooling more gas than GADGET, as is clearly 
evident by the significantly different radial streaming velocities. 
While the lower central temperatures and higher central densities 
in AREPO may suggest a higher cooling rate overall, note that this 
is really only the case in the innermost regions of the haloes, but 
much of the volume and gas mass of the haloes sit at larger radii. 
There, the gas density in AREPO is actually comparable or slightly 
lower than in GADGET, as is also expected due to the reduced re- 
maining gas content we have measured for these haloes in AREPO 
at late times (because more stars have been formed). This means 
that the total energy losses from emitted cooling radiation are actu- 
ally quite comparable for the two codes, especially in the relevant 
region around the cooling radius in the outer parts of the haloes. 
Given these cooling losses, the surprising thing is then actually not 
that AREPO cools out so much gas, it is that GADGET cools out so 
little gas, especially in ^ 10^^ /i^^M© haloes where the cooling 
times are short. 

We argue that this is due to the difference in dissipative heating 
between the two codes that we measure for such quasi- stationary 
haloes in Paper I. In SPH, a combination of viscous dissipation 
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Figure 8. Density, temperature, entropy and radial velocity profiles of~ lO^^/i^^M© haloes at z = 0. We only include gas that is outside the two-phase 
star-forming medium and has temperature T > 100, 000 K, in order to concentrate on properties of hot halo gas. The radius is expressed in units of the virial 
radius (i?20o), temperature is scaled to the virial temperature, density is expressed in units of the mean enclosed baryonic density of our haloes, p200, and 
velocity is scaled to the circular velocity at R200. Red and blue lines are for haloes simulated with AREPO and GADGET, respectively, while different line 
styles indicate different numerical resolutions, as labeled in the density panel. The halo mass range stacked for the measurements is given in the temperature 
panel and includes between 14 and 20 haloes. 



of sonic velocity noise of the gas and efficient damping of turbu- 
lent gas motions injected in the infall regions leads to a significant 
heating of the outer parts of haloes, partially offsetting the cooling 
losses there. In contrast, AREPO produces less dissipative heating 
in this region, so that a stronger cooling flow and higher radial in- 
fall velocities result. The adiabatic heating of the gas on its way 
inward is then maintaining roughly constant or even rising tempera- 
tures for the gas, despite the cooling, until it eventually drops catas- 



trophically in temperature at very small radius, and is then rapidly 
flowing onto the star-forming phase. In addition, stripping of more 
gas-rich disks in AREPO and subsequent mixing of this low-entropy 
material (which is suppressed in GADGET) can further strengthen 
cooling in the central region of AREPO haloes. 
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3.7 Characteristic size and specific angular momentum of 
galaxies 

As an illustration of differences in the sizes and shapes of galaxies, 
we compare in Figure [T0| the gaseous and stellar morphology of a 
typical late-type galaxy formed by our two simulation codes. Each 
projection shows the gas or stellar mass density within a box of 
physical side-length 50 /i~^kpc centered on the galaxy. The galaxy 
shown is one of the most massive galaxies in our simulated volume 
at z = 2, and the same object was selected for a discussion of its 
large-scale environment in Paper I. 

The galaxy resides in a halo of mass Mhaio ^ 9.1 x 



lO^^/i~^M0, and its most important quantitative measurements 
are summarized in Table [2] The visual impression already indi- 
cates that the galaxy is more massive, more extended and has a 
much higher gas fraction in AREPO, which is born out by our quan- 
titative measurements. In fact, the star formation rate is a factor 
of ^ 2 higher in AREPO than in GADGET, illustrating the princi- 
pal difference between the codes we discussed in the previous sec- 
tions. In addition, the gaseous disk is clearly much more extended 
in AREPO, as well as being thinner and more well-defined. While 
the difference in the half-mass radius is ^ 60%, the net specific 
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STARS 





X [kpc/h] 

Figure 10. Example of a massive disk galaxy at redshift z = 2 in AREPO and GADGET L20n512 simulations. The upper two rows of panels show the 
projected surface density of gas in a 50 /i~^kpc physical box centered on the galaxy for two different lines-of-sight. The galaxy is clearly more extended 
in AREPO and has a smoother distribution of dense gas. Some quantitative properties of the object are summarized in Table |2] The edge-on projection in 
AREPO reveals a structure that looks more disky. The lower two rows of panels show projected stellar surface density in the same volume and for the same 
lines-of-sight. Clearly, also the stellar component of this large galaxy has a more disky morphology in AREPO compared with GADGET. Note that gaseous 
and stellar disks in the edge-on projection exhibit higher degree of alignment in the moving-mesh code. 
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Figure 11. Example of a disk galaxy at redshift z = in AREPO and GADGET L20n512 simulations. The upper two rows of panels show the projected 
surface density of gas in a 50 /i~^kpc physical box centered on the galaxy for two different lines-of-sight. The galaxy is clearly more extended in AREPO 
and has a smoother distribution of dense gas. Some quantitative properties of the object are summarized in Table|2] The edge-on projection in AREPO reveals 
a structure that looks more disky. The lower two rows of panels show projected stellar surface density in the same volume and for the same lines-of-sight. The 
stellar component of this large galaxy also has a more disky morphology in AREPO compared with GADGET. 
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Table 2. Properties of the example galaxies shown in Figures[To|and[TT]as measured in our L20n512 simulations. In column (1) we list the simulation name, 
column (2) gives the galaxy baryonic mass, and column (3) the gas mass of our selected galaxy. In column (4) we give the effective half-mass radius of gaseous 
disk in physical units. Column (5) indicates the star formation rate of this galaxy, while columns (6) and (7) give the specific angular momentum of gas and 
stars. Finally, in column (8) we give the angle between the angular momentum vectors of galactic gas and stars. We note that this galaxy has a much higher 
gas fraction, gas extension and a higher star formation rate in AREPO than in GADGET, which is a typical systematic difference we find for massive galaxies 
between the codes. 



T 




Mhcio [Mg /h] M,,„ [Mo /h] 

Figure 12. Dependence of the size of gaseous disks of central galaxies in GADGET and AREPO simulations on the parent halo mass (for a bins size of 0.4 dex 
in halo mass). The solid lines show the average "edge" radius of the dense gaseous disk, i.e. the largest radius within which the azimuthally averaged column 
density of cold dense gas projected onto the disk plane is greater than Nh = 2 x lO^^cm"^. Shaded regions show the 25-75% range of the distribution of 
disk sizes at a given mass, while the dashed lines indicate the standard deviation around the mean. It is apparent that the dense, continuous gaseous disks in 
AREPO are much more extended than their counterparts in GADGET. This holds at all halo masses, both dX z = 2 and z = 0, and it is most dramatic in 
haloes around 10^^ /i^^M© at 2; = 0. 



angular momentum of the gaseous disk is almost 2 times larger in 
AREPO than for the same galaxy in GADGET. 

Similarly, the stellar disk is somewhat more extended and thin- 
ner in AREPO, as seen in Figure [T0| It also has a higher rotational 
support, as indicated by the net stellar specific angular momentum 
which is about 60% higher in AREPO compared to GADGET. This 
is consistent with the visual impression of a more disky stellar dis- 
tribution in AREPO. It is also evident that the stellar and gaseous 
disks are better aligned in AREPO, as measured by the angle be- 
tween the angular momenta of gas and stars in the galaxy. We note 
that our findings in this halo represent typical differences in the 
properties of massive galaxies, and are by no means restricted to a 
single system. 

Figure [TT] shows another example of extended gaseous disk, 
this time of a central galaxy in ^ 2.2 x lO^^/i^^M© halo at 2; = 0. 



The figure demonstrates that gaseous disks are even more extended 
at late times and the relative difference between AREPO and GAD- 
GET is now even larger. The example galaxy contains 3 times more 
gas in its large disk, most of which is non-star forming and repre- 
sents a large extended gaseous reservoir. Galactic interactions be- 
tween such disks can lead to strong, extended tidal features such 
as the ones visible in Figure [6] while compact gaseous disks of 
GADGET galaxies are unlikely to produce them in cosmological 
simulations. Differences in the stellar distribution are visible in the 
bottom panels of Figure [TT] with the galaxy appearing more disky 
in AREPO as a consequence of higher specific angular momentum 
of stars in our moving-mesh simulation (see Table |2]). This galaxy, 
however, does show a large bulge component that is similar in both 
simulations (and is a fairly typical example for this halo mass), sug- 
gesting the need for an efficient feedback mechanism that can re- 
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Figure 13. Specific angular momentum of galactic gas (left panels), galactic stars (middle panel) and for the whole galaxies (right panels), atz = 2 and z = 0, 
as a function of the parent halo mass. The measurements demonstrate that in AREPO (red lines) both the gaseous and stellar components have a higher specific 
angular momentum than in GADGET (blue lines), indicating diskier, more rotationally supported baryonic configurations. Shaded regions and corresponding 
dashed lines show 25-75% of the angular momentum distribution in a 0.4 dex mass bin at 2; = 2, and 0.8 dex mass bin at z = 0. 



move a large fraction of low-angular momentum baryons (e.g. |van| 
|den Bosch|200T||Maller & Dekel|2002[ [Brook et al.|2011| ). 

In order to demonstrate this point more systematically we 
show a measure of the extent of gaseous disks in our highest res- 
olution simulations in Figure [12] For this purpose we project the 
mass of cold, dense gas onto the gaseous disk plane of a central 
galaxy. Next, we find the azimuthally averaged surface density of 
the cold, dense gas in 0.4 kpc wide rings. We define "disk radius" 
as the radius at which the ring- averaged column density drops to 
a column density characteristic of damped Lyman-g absorbers. 
Nh 



12 



it is clear 



2 X 10 cm~ . From the results in Figure 
that galaxies in our high-resolution AREPO simulation are much 
more extended than galaxies in GADGET. Ai z = 2, the difference 
in the average extent of gaseous disks is typically around a factor 
of 1.5 while for the most compact galaxies it is closer to fac- 
tor of ~ 2. In addition, visual inspection indicates that the disks in 
AREPO are not only more extended but are also smoother, thinner 
in edge-on projections and have better aligned gas at different radii 
(as visible in Figure [To]). This makes their appearance more similar 
to the observed galactic disks. In contrast, GADGET galaxies show 
a more irregular gas distribution in the disk, and numerous dense 
gaseous clumps in the galactic surroundings. 

At 2; = 0, below a halo mass of 10^^ /i~^M0, AREPO galaxies 
are also a factor of ^ 1.5 larger, but in more massive haloes the size 
difference is even more dramatic than at 2; = 2. For about an order 
of magnitude in halo mass around 10^^ /i^^M©, galaxies formed 



in our moving-mesh simulations are systematically larger than in 
SPH by a factor of 3-5. Notice that this is a mass range of great 
interest as it corresponds to the halo masses of massive late-type 
disk galaxies such as our own Milky Way. 

The measurement of gas disk sizes, as shown in Figure [12] 
has been specifically chosen to capture as closely as possible visual 
impressions from Figure[TO]and[TT](see also Figures 2 and 4 in Pa- 
per I). The trends of Figure [12] are fully consistent with a more in 
depth analysis presented in Torrey et al. (2011 ), which is based on 
detailed Sersic and exponential disk profile fits. We note that stellar 
half-mass radii are only barely resolved even in our highest resolu- 
tion simulations, being of the order of 1 — 3 times the gravitational 
softening length. We can therefore not make definitive statements 
about the radial extent of the stellar component, based on a rather 
coarse diagnostic such as half-mass radii. Nonetheless, the full ex- 
tent of the stellar component is well resolved in our simulations and 
can be used to infer at least some of the stellar disk properties, such 
as total angular momentum (see below). 

To characterize the degree of orderly rotation of stellar and 
gaseous components of galaxies we calculate the angular momen- 
tum of individual gas cells/particles and stellar particles with re- 
spect to the galaxy center and find the total angular momentum vec- 
tors for gas and stars in a galaxy. We normalize this by the gaseous 
and stellar galaxy mass, respectively, before calculating the mag- 
nitude of the net specific angular momentum. We apply the same 
procedure for the galaxy as a whole. The results are shown in Fig- 
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ure [13] in terms of the net specific angular momenta of gas, stars 
and total baryons (gas+stars) in galaxies in our highest resolution 
simulations. Not too surprisingly, for the cold, dense galactic gas, 
the structural differences we have observed earlier in the character- 
istic radii are also reflected in the specific angular momentum of 
the galactic gas, which is higher in AREPO than in GADGET by a 
factor of ^ 2, explaining the larger sizes of galaxies forming in our 
moving-mesh simulations. 

Substantial differences are also visible in the distribution of 
the stellar specific angular momentum shown in Figure[T3] Because 
galaxies in AREPO are diskier, their specific angular momentum is 
much higher than in GADGET galaxies. This statement applies to 
both z — 2 and z = 0, with typical angular momentum differences 
at 10^^ H-^Mq of the order of 1.5 - 2. 

In Figure [13] we also indicate the spread of the distribution in 
specific angular momentum by indicating the 25% and 75% per- 
centiles of the distribution at a given mass. The distribution of 
specific angular momentum at a given mass is much narrower in 
AREPO than in GADGET. In terms of specific angular momentum 
of both gas and stars in the 10^ — lO^^/i^^M© mass range, there 
is very little overlap between galaxies in AREPO and GADGET. 
For example, galaxies within the lowest 25% of the distribution 
in AREPO often have similar specific angular momentum as the top 
25% of galaxies in GADGET. 

Interestingly, the differences between the net specific angu- 
lar momentum for the total baryonic content of galaxies are even 
larger than for the stellar content, where AREPO galaxies have up 
to a factor of ~ 2 higher specific angular momentum than GAD- 
GET galaxies. This results from a combination of higher gas con- 
tent and higher net specific angular momentum in AREPO, and a 
better alignment of the gaseous and stellar components. 

In our measurements, we see a trend of an increase of the 
gaseous disk extent and the specific angular momentum of galaxies 
with halo mass and time. Both effects are consistent with expecta- 
tions of simple models of the formation of galactic disks ( [Fall &| 
|Efstathiou| 19801 |Mo et al.|1998||Soinerville et al.|2008| ). In fact, 
these trends are qualitatively easily understood in the inside-out 
model of disk formation ( [Fall & Efstathiou|1980| ) where the galaxy 
sizes will approximately scale like Rd oc XRvir, where A is the 
dimensionless spin parameter (> Peebles|[l969| ) and Rd is the disk 
scale length. The specific angular momentum then approximately 
scales as j oc A X i?vir x Vc. Because the spin parameter is approx- 
imately independent of mass and redshift (e.g. Barnes & Efstathiou 
|1987| ), the sizes and specific angular momenta of galactic gas both 
increase with time and mass. Both R^iv and Vc increase with mass, 
causing the dependence of j on halo mass to be stronger than for 
the disk radius, correctly approximating the trends in our simula- 
tions. For a more detailed quantitative model, one however has to 
account for the mass and redshift dependence of the halo concen- 
tration and the angular momentum distribution within haloes. Also, 
one has to postulate conservation of angular momentum of the gas 
during infall, an assumption we do not have to make in our detailed 
hydrodynamical simulations. 

Owing to the increasing specific angular momentum as a func- 
tion of halo and galaxy mass, one may wonder if a large part of 
the differences in Figure [13] is caused by the mass difference of 
galaxies inhabiting halos of the same mass. As we show in Sec- 
tion |3.3| in our highest resolution simulations the galaxy masses 
at a fixed halo mass are very similar between the two codes for 
Mhaio ^ lO^^M©. This means that trends with halo mass also di- 
rectly reflect trends with galaxy mass. However, one needs to be 
more careful at higher masses where there is a systematic offset in 



galaxy masses between AREPO and GADGET at a given halo mass. 
In Figure [14] we show the specific angular momentum for galac- 
tic baryons as a function of galaxy mass. In addition, we show a 
matched sample of galaxies in the mass range 10^° — 2 x 10^^ M©, 
connected with dashed lines to indicate corresponding differences 
in specific angular momentum for galaxies in the same collapsed 
structures. It is clear that the specific angular momentum of AREPO 
galaxies is systematically larger than for GADGET at a given galaxy 
mass. This trend is even clearer for our matched sample: for the 
large majority of galaxies, AREPO galaxies have much higher spe- 
cific angular momentum than galaxies of the same mass in GAD- 
GET. In other words, the baryons that make galaxies in our mov- 
ing mesh simulations either retain more angular momentum or are 
forming from baryons with systematically higher angular momen- 
tum. 

Finally, we note that by looking at measurements of disk sizes 
as a function or resolution, we find that there appears to be a qual- 
itative difference between our simulation techniques with respect 
to the sensitivity of predicted disk sizes on numerical resolution. In 
both codes, the sizes of the galactic gas distributions spuriously in- 
crease with poorer resolution, but this effect is somewhat stronger 
in AREPO. This basic trend is to be expected as poor spatial reso- 
lution will tend to lower the central densities in galaxies and make 
them puffy. In SPH, part of this trend is likely offset by the angu- 
lar momentum transport caused by the artificial viscosity in SPH, 
which will tend to makes galaxies more compact, as we discussed 
in Section[3H 

To summarize, gaseous disks of AREPO galaxies are system- 
atically larger than their counterparts in GADGET. AREPO galaxies 
show a higher degree of rotational support as measured by their 
specific angular momentum, indicating that they are more disky 
and less dispersion-dominated. These differences are particularly 
prominent for the range of halo masses in which large disk galax- 
ies are expected to form in the real Universe. We note that previous 
attempts to model disk galaxy formation in cosmological simula- 
tions with SPH often led to a too low specific angular momentum 
of simulated galaxies (e.g.' Navarro & Steinmetz|200 ). Higher res- 
olution simulations of individual galaxies with SPH hinted that the 
problem might be less severe when the resolution is better (e.g. 
Governato et al.|2004||200 7 ). However, this has not been demon- 
strated yet for a large population of galaxies. Furthermore, we show 
that our moving-mesh simulations provide a dramatic improvement 
over corresponding SPH results and have much higher specific an- 
gular momentum at a given galaxy mass, facilitating formation of 
late- type galaxies. 



4 DISCUSSION AND CONCLUSIONS 

In this study, we have systematically compared results obtained for 
cosmological simulations of galaxy formation with the SPH code 
GADGET against the new moving-mesh code AREPO. The simula- 
tions started from identical initial conditions, employed the same 
high-resolution gravity solver, and used a completely equivalent 
modeling of the physics of radiative cooling, star formation and 
associated feedback processes, the latter being incorporated into 
a subresolution model for the multi-phase structure of the ISM. 
Hence the outcome of our simulations yields a direct test of the 
impact of hydrodynamical solvers on the predictions of cosmolog- 
ical galaxy formation simulations. 

In Paper I of this series, we have focused on a comparison 
of the global baryonic properties predicted by the two simulation 
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Figure 14. Specific angular momentum of galaxies (stars+gas) at z = 2 and z = 0, as a function of the galaxy mass. The measurement demonstrates that in 
AREPO (red lines) galaxies at a given mass have higher specific angular momentum than in GADGET (blue lines). Thick lines show the medians in 0.4 dex 
mass bins at 2; = 2, and 0.6 dex mass bin sat z = 0. Symbols connected with thin dashed lines show a matched sample of galaxies in AREPO (red, diamonds) 
and GADGET (blue, triangles). The matched sample includes the first 20 galaxies from Figures 2 and 4 in Paper I, excluding galaxies undergoing major 
mergers, where the angular momentum changes rapidly. It is clear that for most of the galaxies in our simulations differences in specific angular momenta of 
galaxies are well above the trends expected from their mass differences. 



techniques, while here we consider properties of individual galax- 
ies. Reassuringly, we find that many of their properties are in fact 
in broad agreement, for example the galaxy mass functions are 
quite similar at low masses. However, when scrutinized in detail, 
and especially towards later times and in high mass systems, we 
find rather significant systematic differences in the predicted galaxy 
properties. 

Specifically, massive haloes in AREPO host galaxies that are 
larger, diskier, more massive, and have a higher star formation rate 
than their GADGET counterparts. In the extreme end of the galaxy 
mass function, AREPO can accumulate up to a factor of ^ 2 more 
mass in central galaxies than GADGET. The moving-mesh galaxies 
also have a higher specific angular momentum, and the morphology 
of the gas surrounding the galaxies is clearly different, being much 
less clumpy than in the SPH code. These clumpy condensations of 
gas are largely absent in AREPO; instead, comparatively smooth gas 
disks are formed which tend to align better with the stellar disks. 

There appears to be primarily a single driver for the higher star 
formation rates, lower central entropies and higher gas fractions in 
massive galaxies in AREPO. Multiple lines of evidence show (see 
also Paper I) that the hot gaseous haloes formed in AREPO cool out 
gas more efficiently than their counterparts in GADGET, causing in 
turn much higher star formation rates in the central galaxies of these 
haloes. We find that these different effective cooling rates occur 
despite the fact that the two codes create rather similar gas profiles 
in large haloes at low redshift, which do not suggest an a priori 
cause for a substantially different total cooling emission. In fact, 
we find that the culprit does not lie in a different cooling per se. It 
lies in a different heating. As we have shown in Paper I, the local 
dissipation rates in haloes simulated with SPH and the moving- 
mesh code are quite different. Whereas more entropy (and hence 
dissipative heating) is produced in the infall regions of haloes in 



AREPO compared with GADGET, the opposite is true for most of the 
halo volume inside the virial region, and in particular in the outer 
parts of haloes where the cooling radius for many of these haloes 
lies. This dissipative heating offsets part of the cooling losses in 
SPH, reducing the strength of the cooling flows. 

There are good reasons to believe that this heating effect in 
SPH is in fact spurious. First of all, as we demonstrate explicitly in 
Paper III, the standard formulation of SPH that we use here tends to 
suppress gas stripping out of infalling satellites, as well as any en- 
tropy mixing with background material. As a result, dense gaseous 
clumps penetrate further into haloes in GADGET upon infall, heat- 
ing them through ram-pressure interactions at smaller radii than 
corresponding AREPO simulations, where material is stripped and 
mixed already at much larger radii. The stripping processes them- 
selves, as well as the curved accretion shocks around haloes, in- 
ject subsonic and mildly transonic turbulence into the outer parts of 
haloes. However, as ^Bauer & Springel|p012| ) show, subsonic tur- 
bulent cascades are not represented in a physical way in SPH at the 
achievable resolution (unlike in the moving-mesh code). Instead, 
the turbulence is damped efficiently close to the injection scale. 

In addition, there is significant viscous damping of small-scale 
subsonic noise in SPH. The latter is ultimately sourced by errors in 
SPH's gradient estimate (Abel 2011 ), which trigger jittering mo- 
tions of the gas particles on the scale of the SPH kernel ( |Springel| 
|2010b ) . If the artificial viscosity is increased, the amplitude of this 
sonic noise can be lowered, but the total dissipation through this ef- 
fect stays nearly invariant (Bau er & Springel|2012| . The dynamic 
environment around forming haloes can efficiently regenerate the 
SPH noise, such that its constant damping by the artificial viscosity 
constitutes a spurious heating source. Especially in the slow cool- 
ing regime encountered for large haloes at low redshift, this can 
then affect the cooling flows in haloes. We note that this problem is 
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neither easily cured by higher resolution, lower artificial viscosity, 
nor a larger number of SPH neighbors. What might help is a formu- 
lation of SPH that features a more accurate gradient estimate, which 
seems a prerequisite to eliminate the sub- sonic velocity noise. 

The differences in cooling induce substantial differences in 
the SFRs of galaxies, which are largest for halo masses around 
Mhaio ^ 10^^ /i~^M0 at z = 0, i.e. approximately at the mass 
scale of our own Milky Way galaxy. At this mass scale, we also find 
substantial differences of the predicted angular momentum content 
of both gas and stars, which is reflected in more extended disks 
galaxies forming in our moving-mesh simulations when compared 
to SPH. The differences in the gaseous disk sizes are apparent at all 
halo masses, but around 10^^ /i^^M© they are particularly large 
and reach a factor of 3. Because the angular momentum content 
of halo gas increases with time (see section [3J] l, a more efficient 
gas supply at late times will lead to a relative increase of the spe- 
cific angular momenta and sizes of galaxies. This suggests that the 
efficient late time cooling is not only increasing galaxy masses but 
is also a significant contributor to the differences in galaxy sizes, 
especially at mass ranges where the cooling differences are large. 
However, our results in Figure [14] show that this is not the only 
cause. 

Also, we note that even low-mass galaxies have more ex- 
tended gas distributions and higher gas content in the moving-mesh 
simulations when compared to SPH. These properties together with 
more efficient ram pressure stripping in moving-mesh runs (as we 
show in Paper III), will lead to a much larger mass loss from low- 
mass galaxies as they infall into more massive halos. This means 
that the fraction of the infalling material that will loose angular 
momentum via dynamical friction to dark matter halo will be much 
smaller in AREPO leading to a relative increase of specific angular 
momentum of galactic baryons with respect to GADGET. 

The results we have presented in this paper hence clearly show 
that it does matter for simulations of galaxy formation which hydro- 
dynamical solver is used. We find significant systematic offsets in 
the predicted properties of galaxies formed in GADGET and AREPO 
at all mass scales, but the most striking differences occur in haloes 
around 10^^ /i^^M©. Here the more accurate hydrodynamics 
provided by AREPO leads to late-type galaxies that are large, disky, 
gas-rich, and have high star formation rates. Such galaxies, espe- 
cially at late times, simply appear not to form in cosmological SPH 
simulations when the same physics is included. Our findings also 
suggest that AREPO is able to significantly alleviate the low an- 
gular momentum problem of galaxies formed in cosmological hy- 
drodynamic simulations, which has been a long-standing issue in 
galaxy formation for the past two decades (e.g. |Katz|1992[[Navarro| 
[& Steinmetz^ l997 2000 ). Strong feedback mechanisms, however, 
appear to be needed to both regulate galaxy masses (e.g. |Keres etal.| 
|2009b| ) and remove the low angular momentum baryons that form 
excessive bulges (e.g. "van d en Bosch|[200T1 ). Given the extended 
gaseous disks of galaxies in AREPO and diskier stellar morphol- 
ogy when compared to SPH results, it is likely that the strength 
and type of the required feedback discussed in previous SPH work 
needs to be revisited. This emphasizes that an accurate numerical 
treatment of hydrodynamics is a key factor for successful cosmo- 
logical simulations that aim to form realistic galaxies in a complex 
environment. 
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